%% Code description
% This code runs counterfactual experiments. 
% Figure 6, Figure 7, Table 3, Table 4, Table A.6 and Table A.8 are
% generated.

%% Housekeeping
clear;

% root directory: change accordingly
root_dir = '...';

% add the files in demand & supply estimation folder
addpath(fullfile(root_dir, 'matlab/demand_estimation_AK'));
addpath(fullfile(root_dir, 'matlab/supply_estimation'));

% folder where output figures/tables will be saved
export_tex = 1; 
folder_tex_cf = '...';

%% Load data, demand & supply estimates, and other parameters of the model
load(fullfile(root_dir, 'matlab/supply_estimation/supply_estimates.mat'));

disp(' ')
disp('***************************************************')
disp('Counterfactual starts')
disp('***************************************************')
disp(['alpha= ', num2str(alpha)])

% compute observed insured rate at the market level
D_enroll_st = zeros(N.st,1);
for mm=1:N.st
    ind=find(D.mkt==mm);
    D_enroll_st(mm) = sum(D.s1(ind));
end

% rename some variables
Delta_xi_all = D_unmerge.Delta_xi;
cost1_e_all  = cost1_e; % min loss ratio reg cost coefficient, NNx 1
cost2_1_all  = cost2_1; % rate stability reg cost coefficient, NN x KM
cost2_0      = cost2_0_mat; % rate stability reg fixed cost, NN x KM

clear Delta_xi cost1_e cost2_1 pr_stay

temp = load(fullfile(root_dir, 'matlab/demand_estimation_AK/demand_estimates.mat'));
use_log_n = temp.use_log_n;
clear temp

major = D_unmerge.major;
pr_s = D_unmerge.pr_s;
Ecost_rs = sum(pr_s.*cost_rs,2);
Ecost_rs_pe = Ecost_rs./D_unmerge.s1; % per-enrollee cost

% baseline estimates for rate adjustment cost coeff
Ecost2_1_base = sum(cost2_1_all.*D_unmerge.pr_s,2);

% needed to unmerge prices to provide as initial values
N_st = N.st;
D_major = D.major;
D_nF = D.n_insurers(D_major==0); % N_st x 1, contains # of fringes in a given mkt
arg_share=0;

% in the same order as D, store expanded agg states and rate increases
for jj=1:N.ist
    id_jj = D.id_ist(jj);
    ind = find(D_unmerge.id_ist==id_jj);
    ind = ind(1);
    D.pr_s_m(jj,:) = D_unmerge.pr_s(ind,:);
    D.p2_m(jj,:) = D_unmerge.p2(ind,:);
end

% needed to store in arrays with nrows=N.ist, preserving initial ordering in D
id_ist_old = D.id_ist;

% supply data with expanded agg states + unmerged fringes
D_input = D_unmerge;

% cost scales
scale_vec = [0.5; 0.7; 0.9; 1; 1.3; 1.7; 2]; 

%% CF rate adjustment cost under baseline Medicaid
medi_cf=0; % baseline Medicaid
[ssolve_market_all, ssolve_market, nnon_converge_all, nnon_converge]...
    = cf_rate_adjustmet_cost(medi_cf, scale_vec, cost2_0, cost2_1_all, D_input, D, N,...
    D_major, arg_share, D_nF, N_st,...
    D_enroll_st, ...
    alpha, ...
    Delta_xi_all, cost1_e_all, mu_f_all, sigma_f_all,...
    K, pr_y, N_y, resource_y1, resource_y2, y_medi, x_fc, pr_fc, ...
    rho, crra, sigma_c1, beta_c, beta_f, T1, B1_c, B2_c, B1_f, B2_f, delta, ...
    lr_target, N_f, entry_exitflag_vec, id_ist_old, myopic, Er, cost_weight, use_cost_weight);

%% CF rate adjustment cost under generous Medicaid
medi_cf = 1; % counterfactual Medicaid
y_medi_cf = 2*y_medi; % increase consumption floor

[ssolve_market_all2, ssolve_market2, nnon_converge_all2, nnon_converge2]...
    = cf_rate_adjustmet_cost(medi_cf, scale_vec, cost2_0, cost2_1_all, D_input, D, N,...
    D_major, arg_share, D_nF, N_st,...
    D_enroll_st, ...
    alpha, ...
    Delta_xi_all, cost1_e_all, mu_f_all, sigma_f_all,...
    K, pr_y, N_y, resource_y1, resource_y2, y_medi_cf, x_fc, pr_fc, ...
    rho, crra, sigma_c1, beta_c, beta_f, T1, B1_c, B2_c, B1_f, B2_f, delta, ...
    lr_target, N_f, entry_exitflag_vec, id_ist_old, myopic, Er, cost_weight, use_cost_weight);

% load results under baseline Medicaid
temp = load('temp_variable_cost.mat');
OO = temp.OO;
clear temp

% load results under generous Medicaid
temp = load('temp_variable_cost_medi_cf.mat');
OO2 = temp.OO;
clear temp

% save
save('cf_results.mat')

%% Report results for CF rate adjustment cost under baseline & generoud Medicaid
% Figure 6, Table 4, Figure 7 of the main text are generated. 
shell_report_rate_adjustment(scale_vec, ...
        ssolve_market_all, nnon_converge_all, ssolve_market, nnon_converge,...
        D, K, N, D_enroll_st, pr_y, alpha, B1_c, beta_c, T1, B2_c,...
        ssolve_market_all2, nnon_converge_all2, ssolve_market2, nnon_converge2,...
        lr_target, export_tex, folder_tex_cf, OO, OO2, resource_y1, T2);

close all

%% CF where no states adopt RSR 2000 vs. all states adopt RSR 2000
% get regression coefficients from running 
% supply_estimation/cost_estimates_regression.do in Stata
fc_multiplier = 0.24; % 0.0386/0.1599
vc_multiplier = 1.94; % 0.01357/0.007

OO_rsr2000 = cf_rsr2000(fc_multiplier, vc_multiplier,...
    cost2_0, cost2_1_all, D_input, D, N,...
    D_enroll_st, ...
    alpha, ...
    Delta_xi_all, cost1_e_all, mu_f_all, sigma_f_all,...
    K, pr_y, N_y, resource_y1, resource_y2, y_medi, x_fc, pr_fc, ...
    rho, crra, sigma_c1, beta_c, beta_f, T1, T2, B1_c, B2_c, B1_f, B2_f, delta, ...
    lr_target, N_f, entry_exitflag_vec, id_ist_old, myopic, Er, cost_weight, use_cost_weight, root_dir);

% Report results.
% Table A.8 of Online Appendix is generated
shell_report_rsr2000(OO_rsr2000, ...
     D, K, N, D_enroll_st,  pr_y, alpha, B1_c, beta_c, T1, B2_c,...
     export_tex, folder_tex_cf, fc_multiplier, vc_multiplier, resource_y1, T2)

%% CF commitment: only one price for each contract 
OO_commit = cf_commitment(cost2_0, cost2_1_all, D_input, D, N,...
    D_enroll_st, ...
    alpha, ...
    Delta_xi_all, cost1_e_all, mu_f_all, sigma_f_all,...
    K, pr_y, N_y, resource_y1, resource_y2, y_medi, x_fc, pr_fc, ...
    rho, crra, sigma_c1, beta_c, beta_f, T1, T2, B1_c, B2_c, B1_f, B2_f, delta, ...
    lr_target, N_f, entry_exitflag_vec, id_ist_old, myopic, Er, cost_weight, use_cost_weight, root_dir);

% Report results.
% Table 3 of main text is generated. 
shell_report_commit(OO_commit, ...
     D, K, N, D_enroll_st,  pr_y, alpha, B1_c, beta_c, T1, B2_c,...
     export_tex, folder_tex_cf, resource_y1, T2)

%% Robustness: vary the initial rate regulatory cost
low_multiplier  = 0.5;
high_multiplier = 2;

OO_minloss = cf_minloss(low_multiplier, high_multiplier, ...
    cost2_0, cost2_1_all, D_input, D, N,...
    D_enroll_st, ...
    alpha, ...
    Delta_xi_all, cost1_e_all, mu_f_all, sigma_f_all,...
    K, pr_y, N_y, resource_y1, resource_y2, y_medi, x_fc, pr_fc, ...
    rho, crra, sigma_c1, beta_c, beta_f, T1, T2, B1_c, B2_c, B1_f, B2_f, delta, ...
    lr_target, N_f, entry_exitflag_vec, id_ist_old, myopic, Er, cost_weight, use_cost_weight, root_dir);

% Report results.
% Table A.6 of Online Appendix is generated. 
shell_report_minloss(OO_minloss, ...
     D, K, N, D_enroll_st,  pr_y, alpha, B1_c, beta_c, T1, B2_c,...
     low_multiplier, high_multiplier, export_tex, folder_tex_cf, resource_y1, T2);





